An atomic colour superfluid via three-body loss 
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Large three-body losses in a three-component Fermi gas confined in an optical lattice can prevent 
the occupation of a lattice site by three atoms. This effective constraint not only gives rise to a sup- 
pression of actual three-body loss, but stabilises BCS pairing phases by suppressing the formation of 
trions. We study the effects of the constraint using bosonisation and density matrix renormalisation 
group techniques (DMRG). We discuss the case of lithium experiments, and study the dissipative 
dynamics including loss using time-dependent DMRG with quantum trajectories methods. 



Recent developments in the experimental control of de- 
generate Fermi gases with cold atoms [l[ has paved the 
way for the study of three-component Fermi mixtures 
2]. For attractive two-body interactions, these systems 
offer a chance to observe competition between an atomic 
colour superfluid phase with BCS pairing and a trionic 
phase (with three particles on the same site) 0,110,13. 
A key feature of current experiments, though (e.g., with 
lithium), is the large three-body loss rate observed in 
these mixtures. Here we discuss how these losses, which 
are normally undesirable, can give rise to an effective 
three-body hard-core constraint @, [1| when the gas is 
loaded into an optical lattice Not only does this 
suppress actual loss events, but it also acts to suppress 
trion formation, thus stabilising pairing phases such as 
the atomic colour superfluid. The blockade mechanism 
is related to the continuous quantum Zeno effect, and has 
been demonstrated for two-body interactions in an exper- 
iment with Feshbach molecules 1(J 11 1- The constraint 
has particularly striking consequences in f D, where in the 
absence of such a constraint there is no atomic colour su- 
perfluid (ACS) phase, but rather competition between a 
charge-density wave (CDW) and a phase with symmetric 
(on-site) trions (ST) [ij (see Fig. la). We show that a 
constraint prevents ST formation and produces an ACS 
phase, which competes with a CDW and off-site trions 
(OT) (see Fig. lb). 

A three component Fermi gas in the lowest band of an 
optical lattice is described by the Hamiltonian (H = 1) 



-h.c.)+)^U a mi t o-mi !a - + i, (1) 



where denotes a sum over neighbouring sites, 

c \ a , c i.<? are fermionic operators with a species index a = 
1,2,3, TOi.cr = cj^c^o-, J a are the tunnelling amplitudes, 
and U a are the onsite interaction energy shifts. In the 
following, we will consider balanced densities fh a — h/3 
for total mean number density n, and in typical realisa- 
tions we will have equal tunnelling amplitudes J a = J. 
This model is valid in the limit J a , V a n <C u>, with ui the 
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FIG. 1: Qualitative phase diagram for attractive interac- 
tions U < and equal populations n/3 of each component 
in a three-component ID Fermi gas. These are shown in 
the SU(3) symmetric case (where all pairwise interactions be- 
tween different components are of equal strength), a) without 
and b) with a three-body hard-core constraint arising from 
three-body loss. The unconstrained case is characterised by 
competition between symmetric (on-site) trions (ST) and a 
charge-density wave (CDW). The hard-core constraint sup- 
presses trion formation, stabilising BCS pairing in an atomic 
colour superfluid (ACS), which competes with a CDW and 
off-site trions (OT). 



energy separation between the lowest two Bloch bands. 

Three-body recombination will result in decay into the 
continuum of unbound states, i.e., loss from the optical 
lattice. The decay dynamics can be described by a master 
equation in which loss occurs from a singe site occupied 
by three atoms at a rate 73 Q , 
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where is the system density operator with n atoms, 
U = Cj,iCi,2Cj,3) and the effective Hamiltonian is H e g = 
Hjj — i^t\ti/2. If we begin with an initial state not 
involving three body occupation, then via a mechanism 
analogous to the quantum Zeno effect, a large loss rate 
73 S> J will suppress coherent tunnelling that would pro- 
duce triply occupied sites. For large 73/ J, loss occurs at 
an effective rate that decreases as J 2 /"fz in second-order 
perturbation theory. Then, on a timescale where loss 
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can be neglected, the system dynamics is described by 
the constrained Hamiltonian 

H c = VH V V, V = Y[V J = JJ(1 - ,u, A ,u„ ,,„,,;). 



where V is a projector onto the subspace of states with 
at most two atoms per site. We will initially consider the 
physics of this model, then return to the master equation 
in order to test the assumption that loss is small, and to 
investigate time-dependent preparation of states. 

Below we focus on the case of a ID gas, which allows 
us to make quantitative predictions for realistic experi- 
mental parameters using time-dependent density matrix 
renormalisation group (DMRG) methods [l2| to compute 
both ground states and near-equilibrium time-evolution. 
The ID case differs from higher dimensions due to the ab- 
sence of spontaneous symmetry breaking, but provides a 
striking example of the effect of the hard-core constraint 
in that an ACS phase appears that is not present in the 
absence of the constraint. We first determine a quali- 
tative phase diagram for the constrained model in the 
weakly interacting regime using Tomonaga Luttinger Liq- 
uid (TLL) bosonisation techniques, and then go beyond 
this limit by using time-dependent DMRG methods to 
compute the ground state within the constrained Hamil- 
tonian. We treat both the SU(3) symmetric case where 
all interaction constants are equal (as could be realised, 
e.g., with alkaline earth atoms [HI), and the case where 
interactions are unequal (as are typical in Lithium ex- 
periments at low magnetic fields [2]). At the end we 
combine time-dependent DMRG techniques with quan- 
tum trajectory methods to investigate the production of 
an ACS state by computing time evolution for typical 
experimental parameters under the full master equation. 

Bosonisation formalism for the constrained 
Hamiltonian:- In order to produce a qualitative 
phase diagram, we would like to apply the bosonisation 
formalism to our system. It is therefore necessary to im- 
plement the constraint in a more explicit way. We have 
identified an exact mapping of the constrained fermionic 
Hamiltonian He to an unconstrained fermionic Hamil- 
tonian which automatically respects the constraint, 
at the expense of including higher order interactions. 
We introduce projected operators d\ a = (Ilj^i ^') C L> 
di a = (IX73H T~'j) c i<Ti entirely in terms of which we express 
the Hamiltonian. We verify (i) that the operators di„ 
obey fermionic commutations on the subspace where 
at most two atoms occupy any site, and that (ii) the 
Hamiltonian has vanishing matrix elements in the space 
with occupations greater than two and (iii) acts as 
zero on any state in this latter space. Thus, we arrive 
at a fermionic Hamiltonian with built-in constraint, 
which we analyse with standard bosonisation techniques. 
Here we summarise the results, with the calculations 
presented in more detail in a forthcoming work 



continuum version of (d\ a , dia), from which we can con- 
struct a Hamiltonian by taking the linear combinations 
0c = (fa + fa + fa)/V^> which represents collective 
fluctuations of the total density, and (f> s \ — (fa — fa)/V% 
and 4> S 2 — (fa + fa — 2</>3)/v / 6, which represent the spin 
sectors. If we define TLL parameters K a and conjugate 
momentum fields n a corresponding to each field <f) a 
(a G {c, si, s2}), we then obtain 

a— c,sl,s2 ^ ' 

- cosW^nfa-i] cos[V67r</> S 2] - % cos[v / 87rV>siL 

where v = 2aJsin[7m/3] is the Fermi velocity, a is the 
lattice spacing. The coefficients g ss and g s exhibit non- 
trivial dependence on U/J and m CT , as do K a . 

Correlation functions:- To identify the dominant cor- 
relations and determine a phase diagram, we consider 
the behaviour of correlation functions related to different 
ordering: CDW (for which we compute density-density 
correlations C(x) oc (niUi +x )), ACS (with BCS corre- 
lations P a (x) oc (d| dt o . +1 d i+a . )<T tf i+!B)<T+1 ) ) and OT 



(OT(x) ex (t\J i: 
Each correlation is characterized by a scaling dimension 
that determines its algebraic decay. The decay expo- 
nents for CDW and ACS order are independent of a in 
the SU(3) symmetric case and can be related to TLL pa- 
rameters as: T>cdw = 2K c /3 + K s i/3 + K s2 , T>acs = 
K sl +2/(3K c )+l/(3K s2 ). Correlations for OT order have 
multiple contributions in terms of bosonic operators, with 
two distinct decay exponents in the SU(3) case, T>q T = 
3/(2K c ) + 3K c /2 and V\j T = 3/(2K c ) + K c /6 + 4K s2 /3. 
We can then extract values of the decay exponents as 
a function of U/J and n by expanding K a in the weak 
interaction limit. These can be compared to determine 
a phase diagram, and also benchmarked against values 
extracted from numerical simulations (see Table I). 

Phase diagram from bosonisation (equal interactions) :- 
The ID model Hjj without the hard-core constraint has 
been previously studied for attractive interactions us- 
ing a combination of the TLL formalism and DMRG 
methods 0|. The qualitative phase diagram is depicted 
in Fig. la, and involves competition between symmet- 
ric (on-site) trion (ST) order and charge-density-wave 
(CDW) order. In a wide region near the SU(3) sym- 
metric line U a — U < 0, CDW is dominant for higher 
densities and intermediate interactions. A particular fea- 
ture of this system with U < is that a gap appears in 
the entire spin sector, so that in contrast to a two-species 
Fermi gas, ACS correlations decay exponentially (as do 
spin-density wave correlations). 

In the presence of the constraint this picture changes 
substantially (see Fig. lb with small \U/J\). Outside a 
small region at low densities n < 0.2 and intermediate 
interactions, where off-site trionic correlations compete 
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We introduce three bosonic fields 4>,y(x) related to the 
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FIG. 2: a) ACS-correlations P a (x) with (solid line) and with- 
out (dash-dotted line) the three-body constraint as a func- 
tion of distance i on a 40 site lattice, b) - d) Pa-(x) with 
CDW-correlations C(x) and off-site trion correlations OT(x) 
for different values of U , at a density mi = rfi2 = fhz = 0.2, 
with the constraint. In qualitative agreement with bosoni- 
sation, we observe that (b) for weak coupling CDW clearly 
dominates. As the coupling increases (c, d), ACS becomes 
stronger relative to CDW. For the values presented here, ACS 
dominates off-site trions, with OT(x) decaying exponentially 
for strong coupling. 



with CDW, the spin sector is gapless and all correla- 
tion functions decay algebraically. The ST phase is sup- 
pressed by the constraint, and is replaced by a competi- 
tion between off-site trions (OT), and an ACS phase with 
simultaneous pairing of all three pairs of components. 

Phase diagram from bosonisation (unequal 
interactions):- In the case of unequal interactions, the 
spin and charge sectors are coupled by additional terms 

Hcross = ^Cs2dx<t>cdx4> 2 s +^csldx4 > cdx4 , sl+^ssdx4 , sldx4 , s2, 

with A cs2 = V2(2Ui - U 2 - U 3 )/(6ttv), and 
\csi/V6 = (U 3 - U 2 )/(6nv) = -X SS /VS- When 
the interaction imbalance is strong enough to couple 
different sectors, opening of a gap in the spin sector can 
enhance one-channel pairing with respect to the other, 
and we obtain pairing between the species with the 
strongest attractive interaction. 

Ground state computations with time- dependent 
DMRG (equal interactions) :- In order to quantita- 
tively underpin these results and to go beyond the 
weak-coupling regime, we present calculations based 
on time-dependent DMRG methods 12j. In Fig. 2a, 
we see again the striking comparison between the ACS 
correlations in the ground states for H\j and He, 
which exhibit exponential decay without the hard-core 
constraint, and algebraic decay in the presence of the 
constraint. In Figs. 2b-c we show a comparison of the 
correlations corresponding to the ACS (P a (x)), CDW 
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Vci 
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1.71(0.02)/2.04 1.39(0.02)/1.95 3.3(0.3)/3.00 

1.66(0.02)/2.03 1.40(0.03)/1.96 3.3(0.3)/3.00 

1.60(0.02) 1.40(0.03) 3.1(0.4) 

1.4(0.2) 1.4(0.1) 3.0(0.4) 

1.4(0.4) 1.4(0.3) exp 



TABLE I: Exponents for algebraic decay of correlations, com- 
puted for ground states in a system of 40 lattice sites with 
fher = 0.2 by fitting a power law to the periodic peaks of the 
correlations (see Fig. 2). Errors are given in parentheses, an- 
alytic values in the weak coupling limit are given to the right 
of slashes. We find agreement with the qualitative predictions 
of analytic theory, i.e. for weak coupling we start out with 
the CDW-phase dominant, from where T>cdw remains con- 
stant whilst T>acs and Dot decrease with increasing \U/J\. 
We further observe a transition to exponential decay for the 
off-site trions in the strong coupling limit. 



(C a (x)) and OT (OT(x)), in the ground state of the 
model with a three-body hard-core constraint (He)- 
These are presented for symmetric but varying interac- 
tions from weak to strong coupling, U/J G [—10,-0.3]. 
For increasing interactions, the values of P a {x) become 
stronger compared with C a (x) , so that the ACS appears 
to dominate for stronger interactions. While off-site 
trions still show algebraic decay for U/J > —2, they are 
subdominant to the ACS and CDW, and they decay 
exponentially for U/J < —5 (Figs. Ob-d). This is a 
strong deviation from the weak-coupling bosonisation 
results, and could indicate an instability (e.g., towards 
phase separation) , or the appearance of a gap in the dual 
field of the charge sector. However, with system sizes of 
40 lattice sites, we have not observed any evidence of 
phase separation. 

Extracting the exponents of the algebraic decay of 
the correlation functions further confirms this picture, 
as shown in Table 1. The constrained model sees an en- 
hancement of ACS correlations with decreasing U, while 
CDW correlations generally decay faster as U is low- 
ered. In the weak coupling regime we also generally 
observe good agreement with the perturbative values of 
the exponents from TLL-theory. In the strong coupling 
regime (U/J < —5) for the constrained case, the ex- 
ponents T>acs an d ~Dcdw saturate, with Vacs taking 
a value compatible with the TLL prediction f^ics = 
2/(3iQ + 4/3 in the case where 2/(3K c ) -> 0. This is 
consistent with either an instability (K c — * oo), or n c 
becoming massive. 

Ground state computations with time- dependent 
DMRG (unequal interactions):- In the case of asymmet- 
ric interactions we also observe an ACS pairing, but 
with only the two components paired that exhibit the 
strongest interactions. As an example, we consider the 
case of 6 Li, where in Fig. 3a we plot the Hamiltonian 
parameters as a function of magnetic field strength near 
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500-700 G for a fixed lattice depth. From the pairing 
correlations shown in Fig. 3b, we see that at 520 G, 
the only algebraically decaying correlations are those 
corresponding to components 1 and 3, which have the 
strongest interparticle interaction. 

Time- dependent preparation of states:- Considering 
this example, we now return to the full time-dependent 
dynamics including three-body loss, in order to demon- 
strate a method to produce these ACS states in 6 Li. We 
simulate the many-body master equation on 12-24 lat- 
tice sites by combining time-dependent DMRG methods 
with quantum trajectories techniques, as described in Q. 
We assume that the lattice is initially loaded at a mag- 
netic field of 615 G, where the repulsive interactions (see 
Fig. 3a) will stabilise the system in the presence of loss. 
We then consider a time-dependent ramp of the mag- 
netic field to 500 G. The characteristics of the ramp we 
choose (shown in Fig. 3c) are: (i) it is adiabatic until 
565 G, where the components 2 and 3 become paired, 
(ii) it is fast from 565 — 500 G, where onsite trions be- 
come energetically favoured U a < U2) and where 
for fields larger than 520 G, 73 is too small to prevent 
triple occupation (see Fig. 3a), and (iii) after a hold time 
Th,\ — 16J -1 we add a swap between species 1 and 2 
via a fast laser pulse at the end of the ramp, after which 
there is a second hold time T h>2 = 20J -1 . The ACS 
correlations in the final state after the swap (Fig. 3d) 
then exhibit dominant pairing between species 1 and 3, 
as would be expected in the ground state for the param- 
eters of 6 Li at B = 500 G with the constraint. At the 
same time, the other two pairing channels are clearly sub- 
dominant as evidenced by their much weaker correlators, 
while trion formation is also strongly suppressed during 
the ramp. Further, the probability that not even a single 
decay event occurs during the ramp shown here is ~ 80% 
for 12 lattice sites. 

These ideas can be extended to higher dimensions, 
where the ACS should also be stabilised due to sup- 
pressed trion formation. 
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